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Abstract 

The light intensity distribution in a shower image and its implications to the primary 
energy reconstructed by the fluorescence technique are studied. Based on detailed 
CORSIKA energy deposit simulations, a universal analytical formula is derived 
for the lateral distribution of light in the shower image and a correction factor is 
obtained to account for the fraction of shower light falling into outlying pixels in the 
detector. The expected light profiles and the corresponding correction of the primary 
shower energy are illustrated for several typical event geometries. This correction 
of the shower energy can exceed 10%, depending on shower geometry. 



1 Introduction 



One of the methods of extensive air shower (EAS) detection is recording fluo- 
rescence light emitted by nitrogen molecules in the air along the shower path. 
For very high energies of the primary particle, enough fluorescence light is 
produced so that the shower can be recorded from a distance of many kilo- 
meters by an appropriate optical detector system [1,2,3]. As the amount of 
fluorescence light is closely correlated to the ionization energy deposit in air, 
it provides a calorimetric measure of the primary energy. 
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The field of view of a fluorescence detector (FD) telescope is divided into 
many pixels. For example, in case of the Pierre Auger Observatory (PAO) [4] 
each pixel views 1.5° of the sky and records the received light in 100 ns time 
intervals. A shower passing through the telescope field of view triggers some 
pixels, which form together a "shower track". The lateral width of this track 
depends on shower geometry but can well be larger than the pixel size. 

For a precise energy determination one needs to collect the available signal as 
completely as possible, i.e. from all detector pixels which receive light from the 
shower. On the other hand, adding signals from many pixels implies adding the 
background noise as well. Therefore it is important to include in the analysis 
only a small number of pixels which contain the "true" shower signal. 

In this paper, Monte Carlo simulations of the shower image are presented. 
Based on the spatial energy deposit of shower particles as calculated by COR- 
SIKA [5,6] it is shown that the lateral shower spread can be well parameterized 
as a function of the shower age parameter only. The derived parameterization 
can be used for reconstruction of shower profiles. This is illustrated by ap- 
plying the new parametrization to the reconstruction of the primary shower 
energy for several simulated events in a fluorescence detector. 

The plan of the paper is the following: the definition of the shower width and 
algorithm of fluorescence light production based on the CORSIKA simulation 
of energy deposit density are described in Section 2. In Section 3 an analytical 
parametrization is derived and its implementation in energy reconstruction 
procedure is discussed. Conclusions are given in Section 4. 



2 Properties of shower image 

2. 1 Shower width and shape function 

Photons which constitute an instantaneous image of the shower originate from 
a range of shower development stages [7], namely from the surface S shown 
in Figure 1. These simultaneous photons are defined as those which arrive at 
the FD during a short time window At. During this At (corresponding to a 
small change of the shower position in the sky by Ax) the shower front moves 
downward along the shower axis by a small distance Al ~ RA\/ sin(^), where 
R is the distance from FD to the volume AV and Qi is the angle between the 
shower axis and the direction towards FD. This means that the small element 
of surface S corresponds to a small volume AV. The number of photons which 
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dNL ec \ AV = f(X,r)dS ± x x / W(\)r)(\)d\, ( 1 ) 



arrive to the FD from volume AV can be calculated as: 

A 
AttR 2 

where f(X, r) is the distribution of light emitted, dS± is the projection of the 
surface dS onto a surface perpendicular to direction of the shower axis, A is 
the light collecting area of the detector, W(\) is the light transmission factor, 
r](\) is the normalized fluorescence wavelength spectrum. These photons form 
an instantaneous image of the shower which can be described by the angular 
distribution of light recorded by the FD 
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where a is the small angle between the direction to the center of the image 
spot and the direction to volume AV, is the azimuth angle. 

The size of shower image is defined as the minimum angular diameter 2a of 
the image spot containing a certain fraction F(a) of the total light recorded 
by the FD. 

A shower viewed from a large distance has, to a very good approximation, a 
circular image, independent of the direction of shower axis [2,7]. The intensity 
distribution of light in this image is proportional to the lateral distribution 
of the emitted fluorescence light in the shower at the viewed stage of evolu- 
tion. Therefore the fraction of light received F(a) can be obtained from the 
corresponding fraction F(r) of light emitted around the shower axis 

a r 

F(a) = j / 7 (a')27raW ~ J f(r')2irr'dr = F(r), (3) 
o o 



where f(r) is the (normalized) lateral distribution of fluorescence light emit- 
ted. Here we have neglected the fact that photons from the side of the shower 
front facing the detector have been emitted at a time later than photons from 
the farther side of the shower. As will be shown later, relevant lateral distances 
are in the range 150 to 300 m and hence this is also the distance scale for the 
maximum time differences of emission [2]. As there is no significant change 
of the lateral distribution expected of a shower traversing ~300 m of air, the 
effect of different emission is negligible in our calculation. 

In the following we shall consider showers close to the FD. For these showers, 
the optical image size is mainly determined by the geometric size of the shower 
disk. Light absorption and multiple scattering cause only a minor, negligible 
modification of the shower image. 
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The main task is therefore to derive f(r), which is also referred to as the shape 
function, since the brightness distribution of the shower image depends on the 
shape of f(r). At the first approximation, the f(r) is proportional to the num- 
ber of particles in the shower at a given lateral distance, assuming a constant 
fluorescence yield per particle in the shower. For an electromagnetic shower 
this number of charged particles is given by the Nishimura-Kamata-Greisen 
(NKG) function [8]. In this case the F(r) can be determined analytically, as 
has been shown in Ref. [9]. However, F(r) derived in this way does not de- 
scribe the light distribution well in case of hadronic showers. This is due to 
the fact that the number of particles in a hadronic shower does not follow the 
NKG distribution well. 



2.2 CORSIKA approach based on energy deposit density 

Although the assumption that the amount of emitted light is proportional to 
the number of particles is adequate for a determination of the fluorescence 
signal in many cases, it is here not suited for the following reasons. The fluo- 
rescence yield is proportional to the ionization energy deposited by the shower 
rather than to the total number of charged particles. Furthermore, the simu- 
lated number of particles in a Monte Carlo calculation depends on the thresh- 
old energies chosen by the user [10,6], above which particles are simulated. 
Particles falling below the threshold energy are discarded. 

A better approximation for the fluorescence yield can be obtained by using 
the energy deposit dE(X)/dX as a function of atmospheric slant depth inter- 
val dX together with a density- and temperature-dependent fluorescence yield 
Y(p,T) [11,12,13]. In this approximation the distribution of photons emitted 
around the shower axis is proportional to the lateral distribution of energy de- 
posit, f(r) ~ dE }x' r ^ ' wriere dX v = dX cos(#) is the vertical depth interval and 
9 is the shower zenith angle. The distribution of energy deposit dE(X)/dX v 
is calculated with the CORSIKA shower simulation program as the sum of 
the energy released by charged particles with energies above the simulation 
threshold plus the releasable energy fraction of particles discarded due to the 
energy cut. More specifically, the following approximation is used [6]: 

dEi^X^ Ei on i z E e ±^ cu i Ey^ cu t 1 E^±^ cu i ^-^had,cut i j>\ 

dX v ~ AX V AX V AX V 3 AX V 3 AX V U 

where E ioniz is the ionization energy deposit of all charged particles traversing 
the depth interval AX V . E itCUt denotes the energy of particles of type i falling 
below the simulation threshold within this interval. 

In the following we will study the lateral distribution of energy deposit den- 
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sity in air showers, as it is directly proportional to the number of expected 
fluorescence photons. Using CORSIKA a two-dimensional energy deposit dis- 
tribution around the shower axis is stored in histograms during the simulation 
process for 20 different vertical atmospheric depths. Each of the 20 horizontal 
layers has a thickness of AX V = 1 g/ cm 2 and corresponds to a certain atmo- 
spheric depth: the first one to X\ = 120 g/cm 2 and the last one to X20 = 870 
g/cm 2 . Linear interpolation between the observation levels is performed in or- 
der to get the lateral distribution at a given vertical depth X n located between 
two CORSIKA observation levels X k and X k+1 . The fraction of energy deposit 
F(r) is calculated by numerically integrating the histograms up to the lateral 
distance r. 

The shower simulations are performed with the hadronic interaction models 
GHEISHA [17] (for interactions below 80 GeV) and QGSJET 01 [18]. Electro- 
magnetic interactions are treated by a customized version of the EGS4 [19,20] 
code. To reduce computing time, a thinning algorithm [21] is selected within 
CORSIKA. The thinning level of 10~ 6 has been chosen with the so-called op- 
timum weight limitation [22,23]. This ensures that the artificial fluctuations in 
the longitudinal shower profiles introduced by the thinning method are suffi- 
ciently small for this analysis [10]. The kinetic energy thresholds for explicitly 
tracking particles were set to: 100, 100, 0.25, 0.25 MeV for hadrons, muons, 
electrons and photons, respectively. 



3 Results 

The knowledge of the F(r) function can be used to calculate the "true" signal 
(light) from shower, which may be divided among several neighboring detec- 
tor pixels. Below we propose a universal parameterization of F(r) based on 
CORSIKA simulations and show how the "true" signal can be estimated with 
this parameterization. 

3.1 Fraction of energy deposit density from CORSIKA 

In the following we study the dependence of the lateral energy deposit density 
on various variables. 

A natural transverse scale length in air showers, which proves to be useful for 
obtaining a universal parameterization of the lateral distribution, is given by 
the Moliere radius [15] 

r M = E s — , (5) 
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where E s ~ 21 MeV is the scale energy, eo = 81 MeV the critical energy and 
Xi = 37 g/cm 2 the radiation length in air. The local Moliere radius at a given 
atmospheric depth of shower development (at altitude h) can be obtained 
by dividing Eq. (5) by the air density, p(h), and is approximately given by 
r M = 9.6gcnT 2 /p(/i). 

It is also well known that the distribution of particles in a shower at a given 
depth depends on the history of the changes of r M along the shower path 
rather than on the local tm value at this depth. To take this into account, 
the tm value is calculated at 2 cascade units (radiation length Xi) above the 
considered depth [15]. Using the value of the Moliere radius calculated based 
on the atmospheric profile (the US Standard Atmosphere [16]) for vertical 
depth X n — 2Xi cos(#), the fraction of energy deposit density F(r*) versus the 
distance in Moliere units r* = t/tm is found. The knowledge of F(r*) gives 
a possibility to study variation of the shape of energy deposit density due to 
properties of the atmosphere. The variation of the density of the atmosphere 
along the path of a shower affects the Moliere radius and consequently also 
the radial particle distribution. 

To characterize the development stage of a shower, we introduce the shower 
age parameter 

1 + 2X max / 'X' 



where X max is the atmospheric depth of shower maximum extracted from 
simulated data 1 . With this definition, a shower reaches its maximum at s = 1. 

Figure 2 presents the integrals of energy deposit density F(r) and F(r*) for a 
vertical proton shower with primary energy E = 10 EeV, obtained at different 
atmospheric depths. It is seen in Figure 2A that the shape of this integral 
distribution varies considerably only at depths smaller than 360 g/cm 2 . At 
larger depths, and in particular around the shower maximum, the variation 
of integral energy deposit profile is not significant. However, this variation is 
larger when one plots this integral versus distance measured in Moliere units, 
as shown in Figure 2B. 

Figure 3A shows the dependence of the integral of energy deposit density F(r) 
on energy and primary particle. It is seen that the integral profile only slightly 
depends on energy and primary particle. The differences are even smaller if we 
plot the fraction of energy deposit density versus distance in Moliere units, as 
shown in Figure 3B. The same shape of the F(r*) profile for different primaries 

1 X max was determined by fitting a Gaisser-Hillas type function [14] to the COR- 
SIKA longitudinal profile of energy deposit. 



6 



and energies means that variations of the F(r) profile are mainly due to the 
atmospheric effect i.e. dependence of the Moliere radius on altitude. For the 
same shower geometry, there are different altitudes of the maxima of proton 
and iron showers, and in consequence different values of Moliere radius r M . 
Since tm determines the lateral spread of particles in the shower, the shape 
function F(r) becomes broader for iron showers (higher altitude of shower 
maximum than for proton shower). 

Figure 4 presents the dependence of the integral of energy deposit profile 
on zenith angle. We note that CORSIKA energy deposit lateral profiles are 
obtained for horizontal planes at the given observation level, so if one compares 
energy deposit densities between vertical and inclined showers, a projection 
of densities from horizontal plane to the plane normal to the shower axis is 
performed for inclined showers. The corrected profile for a shower inclined at 
9 = 45° is shown in Figure 4A by the solid line. It is seen that this profile differs 
from the profile obtained for a vertical shower. This means that the shape of 
F(r) depends on the zenith angle. This dependence can be understood if 
one takes into account the influence of the atmospheric effect on the energy 
deposit profile F(r). For a homogeneous atmosphere, the shape function for 
inclined and vertical showers must be the same for the same development 
stage because the Moliere radius does not change with altitude. In case of an 
inhomogeneous atmosphere, differences of the shape function between vertical 
and inclined shower should be proportional to the differences of the Moliere 
radius (i.e density of air). Thus if changes of the shape function are due only 
to the atmospheric effect, then F(r*) profile should be the same for vertical 
and inclined showers. Figure 4B confirms this assumption. 

The analysis of Figs. 2B, 3B and 4B leads to the following conclusion: the 
lateral shape of the energy deposit density versus distance from shower axis 
measured in Moliere units is independent of the primary energy, primary par- 
ticle type and zenith angle. It depends, to a good approximation, only on the 
shower age. Figure 5 confirms this conclusion, too. In this Figure we present 
the integral of the energy deposit density for different age parameters for 10 
individual proton and 5 individual iron showers with different zenith angles 
(9 = 0°,45°,60°) and energy 10 EeV. It is seen that the shower-to-shower 
fluctuations are strongly reduced for a given age when we correct F(r) profiles 
for the atmospheric effect i.e. plot F(r*). Also, there are no differences in the 
shape of F(r*) for showers with different zenith angles and primary particle 
type. This means that it is possible to find a universal function which describes 
the shape of energy deposit density as a function of shower age only. Following 
our earlier work [9] we will use the function 

F(r*) = 1-(1 + a(s)r*y bis \ (7) 
where the parameters a(s) and b(s) are assumed to be functions of shower 
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age. Fits of this functional form to the integral of energy deposit density were 
performed for the data from Figures 5B, D, F and are shown in Figure 6. The 
values of the parameters a(s) and b(s) for different shower ages are presented 
in Figure 7. The age dependence of a(s) and b(s) parameters is well described 
by 

a(s) = 5.151s 4 - 28.925s 3 + 60.056s 2 - 56.718s + 22.331, (8) 
b(s) = -1.039s 2 + 2.251s + 0.676. (9) 

Thus, Eqs. (7), (8) and (9) give us a universal function which describes the frac- 
tion of energy deposit density within a specified distance from the shower axis 
for different energies, zenith angles and primary particles. Moreover, Eq. (7) 
can be used to simulate the size of the shower image not only at shower max- 
imum like in Ref. [9], but also for any shower development stage. Inverting 
Eq. (7) and taking into account the distance from the detector to the shower 
(Ro) we can find the angular size of the image a that corresponds to a certain 
fraction of the total fluorescence light signal: 

a(s) = 2arctan (-?-) = 2arctan f -^TT^ 1 ~ F(r))- 1/h{s) - 1) 

3.2 Application in energy reconstruction procedure 

As the shape of the lateral distribution of energy deposit can be well described 
by Eq. (7), it may be used to take into account the knowledge on the shower 
width in the procedure of shower reconstruction in the fluorescence detector. 

One of the first steps in shower energy reconstruction is the calculation of 
the light profile at the aperture of the detector, based on the signal recorded 
by the detector pixels. This signal is converted to the number of equivalent 
photons at the detector diaphragm. For example, one procedure to determine 
such a profile is described in [24]. This algorithm uses as input the recon- 
structed geometry to locate the shower image on the FD telescope camera in 
100 ns intervals. Next, the signal (charge) and noise from pixels lying within 
a predetermined angular distance Xs/n from the instantaneous position of the 
image spot center are collected to find the radius x™/n that maximizes the 
signal-to-noise ratio over the whole shower track. Finally, the charge in each 
100 ns time interval (time slot), L s /N(t), within that radius x's/n ^ s f° un d and 
converted to the number of photons using calibration constants. 

This procedure works well for distant showers, when the light collected within 
the radius Xs/n corresponds to about 100% of the true signal, but some differ- 
ences between the signal within Xs/n an d the true signal may exist for nearby 
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showers. In the following we investigate this problem and estimate a correction 
to the described reconstruction algorithm. 

The necessary shower reconstructions were performed using Flores-Eye [25] 
and FDSim [26] programs. First, FDSim was used to generate events based on 
the Gaisser-Hillas parameterization. Next, the geometrical and energy recon- 
struction were performed using Flores-Eye. The reconstructed geometries for 
some events are listed in Table 1. Using these geometries, we find the collected 
signal L s / N (t) within angular distance Xs'/n &t eacn time interval 2 . Then, for 
the given Xs'/n the effective radius around shower axis r = Ro tan(x^y^) and 
the fraction of light F(r ) based on the function F(r) was calculated. The frac- 
tion F(ro) for events listed in Table 1 is shown in Figure 8. It is seen that for 
Event 1 F(ro) changes from 89% for a distance-to-shower Rq = 7.0 km to 87% 
for R = 6.0 km. For other events, the collected fraction of the signal within 
Xs'/n increases with increasing distance-to-shower and equals on average about 
91%, 94% and 99% for Event2, Event3, Event4, respectively. In other words, 
some portion of the signal, which falls beyond Xs'/n * s missing i n the shower 
reconstruction procedure. To take into account this lost portion of the signal, 
the signal L s / N (t) is rescaled according to formula L true (t) = L s / N (t)/F(r ). 
In this way, one takes into account the shape of lateral distribution of energy 
deposit and obtains the new integrated charge L true (t) for each time slot. Thus 
the part of the signal which was contained in neighboring pixels outside Xs'/n 
is accounted for. 

It should be pointed out that in general, any reconstruction procedure has to 
take into account the pixellation of the detector. Independent of the specific 
approach, the correction procedure developed in this work can be applied. In 
the following, we demonstrate the influence of this correction on the light pro- 
file and on energy determination. In Figure 9A it is shown that our correction 
leads to considerable differences between the L s / N {t) (dashed line) and L true (t) 
profile (solid line) for a nearby shower (Event 1). In case of distant showers 
(like Event4) the profile is almost unchanged (see Figure 9B). It should also 
be noted that changes in the detector-to-shower distance are accounted for in 
this approach. A "differential" correction is applied (i.e. for each time slot), 
which also leads to a better reconstruction of the longitudinal shape (and thus 
X max ) of the shower. 

Accepting only a fraction of the signal contained within Xs'/n directly influ- 
ences the reconstructed primary energy of the shower. In Table 2 we present 
the influence of the L s / N (t) correction on the Gaisser-Hillas fit to the re- 
constructed number of particles in the showers. It is seen that this correction 
changes both the number of particles at the shower maximum and the position 



2 The value of X^/n depends on geometry, but for events listed in Table 1 it equals 
about 1.2°. 
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of the shower maximum. These changes lead to different estimates of primary 
energy. In the last column the relative differences = (-E^™ 6 — E^ N )/E^ 
are listed. One sees that ks is always positive and decreases from 14% for a 
distance to shower maximum of i? =6.5 km to 2% for i? =23 km. 



4 Conclusions 

In this work, the distribution of light in the shower optical image is analyzed, 
based on the lateral distribution of energy deposited by the shower, as derived 
from CORSIKA simulations. The lateral distribution of energy deposited is 
parameterized with a functional form inspired by the NKG distribution. The 
angular distribution of photons arriving simultaneously at the detector (i.e. 
the intensity distribution of light in the instantaneous image of the shower) is 
obtained. The shape of this distribution can be approximated by a universal 
function that depends on the shower age only. 

This universal function is used to derive a correction to the shower energy 
due to the fraction of light falling into detector pixels located far from the 
center of the shower image. In the usual procedure of shower reconstruction, 
signal-to-noise ratio is optimized, so that pixels lying far from the center of the 
shower image are not included in the analysis. The percentage of shower signal 
in those outlying pixels was determined in this paper based on the lateral dis- 
tribution of light in the shower image. The signal recorded by the fluorescence 
detector in the accepted central pixels is rescaled, so that a corrected light pro- 
file of the shower is obtained. For events examined, this correction increases 
the estimated shower energy by 2-14%, depending on the detector-to-shower 
distance. 
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Table 1 

Characteristics of events used for the comparisons in this paper. The shower zenith 
angle 9, azimuth angle cf>, core position X c , Y c are measured relative to FD detector. 



Event 


e 


<t> 


x c 


Y c 




(deg) 


(deg) 


(km) 


(km) 


Event 1 


0. 


-15 


-0.820 


5.860 


Event2 


60 


175 


4.535 


6.906 


Event3 


45 


20 


2.282 


7.066 


Event4 


45 


50 


-7.665 


19.11 



Table 2 

Comparison of Gaisser-Hillas function parameters based on the L true (t) and the 
L S i N {t) light profiles and their influence on primary energy. The number of particles 
at shower maximum N max , corresponding slant depth X max , estimated energy E$ 
and relative difference fee = (£'g r " e — E^ N ) / E^ N are listed as a function of distance 
Rq from detector to shower maximum. Calculations are made using the function 
F(r*) described by Eq. (7). 



Event 


Ro 


N S/N 
1 ^max 


ATtrue 
max 


V S/N 
^■max 


vtrue 
max 


E S/N 


E true 


k E 




(km) 


(io 9 ) 


(10 9 ) 


(g/cm 2 ) 


(g/cm 2 ) 


(EeV) 


(EeV) 


(%) 


Event 1 


6.4 


0.93 


1.06 


701 


706 


1.370 


1.562 


14 


Event2 


8 


6.57 


6.88 


759 


767 


9.853 


10.40 


6 


Event3 


11 


2.12 


2.19 


637 


642 


2.950 


3.100 


5 


Event4 


23 


12.85 


13.10 


752 


753 


19.20 


19.57 


2 
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Fig. 1. Geometry of an EAS as seen by the fluorescence detector. Photons which 
arrive simultaneously at the FD originate from the surface S. See text for more 
details. 
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Vertical proton shower with 10 EeV 
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Distance from shower axis in Moliere units 

Fig. 2. (A) Integral of energy deposit density F(r) versus distance from shower axis 
for proton shower; (B) Integral of energy deposit density F{r*) versus distance from 
shower axis measured in the Moliere units. Different lines described in the top panel 
correspond to profiles obtained for different vertical atmospheric depths. The profile 
corresponding to shower maximum is marked by dots. 
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Fig. 3. (A) Integral of energy deposit density versus distance from shower axis; 
(B) The integral profiles versus distance measured in Moliere units; The profiles are 
shown for vertical showers (at s=l) with different primary particle type and energy. 
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Distance from shower axis in Moliere units 

Fig. 4. (A) Integral of energy deposit density versus distance from shower axis for 
vertical and inclined (9 = 45° ) proton showers; (B) The integral profile measured in 
Moliere units. The profiles are shown for 10 EeV showers at s = 1. 
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Fig. 5. (A), (C), (E) Integral of energy deposit density versus distance from shower 
axis; (B), (D), (F) Integral of energy deposit density versus distance measured in 
Moliere units. Individual profiles are obtained for showers with different age, zenith 
angle and primary particles. 



17 




18 



CO 



-Q 




0.6 0.8 
Shower age s 




0.6 0.8 
Shower age s 

Fig. 7. Values of parameters a(s) and b(s) of Eqs. (8) and (9) obtained based on 
integral of CORSIKA energy deposit density for vertical showers at energy 10 EeV. 
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Fig. 8. Fraction of light collected within the angle Xsyyv versus time for some simu- 
lated events. See text for more details. 
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Fig. 9. (A) Comparison of light profiles versus time slot number for Eventl. The 
Lgij^it) light profile (dashed line) is obtained based on the signal-to-noise algorithm 
implemented in FloresEye. The Lt rue {t) light profile (solid line) is obtained using 
the shape of F(r) function derived from Eq. (7). (B) Comparison of light profile 
versus time slot number for Events . 
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